home *** CD-ROM | disk | FTP | other *** search
/ Amiga Plus 1995 #5 & #6 / Amiga Plus CD - 1995 - No. 5 and 6.iso / pd / serien / purity / nr.49 / astronomie / units / moolib.p < prev    next >
Text File  |  1995-06-24  |  20KB  |  473 lines

  1. {---------------------------------------------------------------------------}
  2. { Unit MOOLIB: Mondbahn                                                     }
  3. {---------------------------------------------------------------------------}
  4.  
  5. UNIT MooLib;
  6.  
  7. {$Projekt Astronomie}
  8.  
  9. INTERFACE;
  10.  
  11.  FROM Astronomie USES MatLib,PnuLib,SphLib;
  12.  
  13. {---------------------------------------------------------------------------}
  14. { MINI_MOON: Mondkoordinaten geringer Genauigkeit (ca.5'/1')                }
  15. {            T  : Zeit in jul.Jahrh. seit J2000  ( T=(JD-2451545)/36525 )   }
  16. {            RA : Rektaszension (in h)                                      }
  17. {            DEC: Deklination (in Grad)                                     }
  18. {            (Aequinoktium des Datums)                                      }
  19. {---------------------------------------------------------------------------}
  20.   Procedure MINI_MOON(T:Real; Var RA,DEC:Real);
  21.  
  22. {---------------------------------------------------------------------------}
  23. {                                                                           }
  24. { MOON: analytische Mondtheorie nach E.W.Brown (Improved Lunar Ephemeris)   }
  25. {       mit einer Genauigkeit von ca. 1"                                    }
  26. {                                                                           }
  27. {       T: Zeit in julianischen Jahrhunderten seit J2000 (Ephemeridenzeit)  }
  28. {          (T=(JD-2451545.0)/36525.0)                                       }
  29. {       LAMBDA: geozentrische ekliptikale Laenge (Aequinoktium des Datums)  }
  30. {       BETA:   geozentrische ekliptikale Breite (Aequinoktium des Datums)  }
  31. {       R:      geozentrische Entfernung (in Erdradien)                     }
  32. {                                                                           }
  33. {---------------------------------------------------------------------------}
  34.   Procedure MOON(T:Real; Var LAMBDA,BETA,R: Real);
  35.  
  36. {---------------------------------------------------------------------------}
  37. { MOONEQU: aequatoriale Mondkoordinaten                                     }
  38. {          (Rektaszension RA und Deklination DEC in Grad, R in Erdradien)   }
  39. {          T in julian.Jahrhndt. seit J2000 ( T:= (JD - 2451545.0)/36525 )  }
  40. {          Die Koord. beziehen sich auf das wahre Aequinoktium des Datums.  }
  41. {---------------------------------------------------------------------------}
  42.   Procedure MOONEQU(T:Real; Var RA,DEC,R: Real);
  43.  
  44. {---------------------------------------------------------------------------}
  45. { T_FIT_MOON: Berechnet die Tschebyscheff-Entwicklung der                   }
  46. {             Koordinaten des Mondes (Reihen fuer RA,DEC und Radius).       }
  47. {                                                                           }
  48. {       TA      : Beginn des Entwicklungsintervalls (jul.Jahrh. seit J2000) }
  49. {       TB      : Ende des Entwicklungsintervalls ( TB < TA + 1 Monat )     }
  50. {       N       : Ordnung der Entwicklung                                   }
  51. {       RA_POLY,DE_POLY,R_POLY: Tschebyscheff Polynome fuer RA,DEC,R        }
  52. {---------------------------------------------------------------------------}
  53.   Procedure T_FIT_MOON(TA,TB:Real; N:Integer;
  54.                        Var RA_POLY,DE_POLY,R_POLY:TPOLYNOM);
  55.  
  56. IMPLEMENTATION
  57.  
  58. {---------------------------------------------------------------------------}
  59. { MINI_MOON: Mondkoordinaten geringer Genauigkeit (ca.5'/1')                }
  60. {            T  : Zeit in jul.Jahrh. seit J2000  ( T=(JD-2451545)/36525 )   }
  61. {            RA : Rektaszension (in h)                                      }
  62. {            DEC: Deklination (in Grad)                                     }
  63. {            (Aequinoktium des Datums)                                      }
  64. {---------------------------------------------------------------------------}
  65. Procedure MINI_MOON;
  66. Const P2     = 6.283185307;
  67.       ARC    = 206264.8062;
  68.       COSEPS = 0.91748;
  69.       SINEPS = 0.39778; { cos/sin(Ekliptikschiefe) }
  70. Var   L0,L,LS,F,D,H,S,N,DL,CB    : Real;
  71.       L_MOON,B_MOON,V,W,X,Y,Z,RHO: Real;
  72. Begin
  73.   { mittlere Elemente der Mondbahn }
  74.   L0:=   frac(0.606433+1336.855225*T); { mittl. Laenge des Mondes (in r) }
  75.   L :=P2*frac(0.374897+1325.552410*T); { mittl. Anomalie des Mondes }
  76.   LS:=P2*frac(0.993133+  99.997361*T); { mittl. Anomalie Sonne }
  77.   D :=P2*frac(0.827361+1236.853086*T); { Diff. Laenge Mond-Sonne }
  78.   F :=P2*frac(0.259086+1342.227825*T); { Knotenabstand }
  79.   DL := +22640*sin(L) - 4586*sin(L-2*D) + 2370*sin(2*D) +  769*sin(2*L)
  80.         -668*sin(LS) - 412*sin(2*F) - 212*sin(2*L-2*D) - 206*sin(L+LS-2*D)
  81.         +192*sin(L+2*D) - 165*sin(LS-2*D) - 125*sin(D) - 110*sin(L+LS)
  82.         +148*sin(L-LS) - 55*sin(2*F-2*D);
  83.   S := F + (DL+412*sin(2*F)+541*sin(LS)) / ARC;
  84.   H := F-2*D;
  85.   N := -526*sin(H) + 44*sin(L+H) - 31*sin(-L+H) - 23*sin(LS+H)
  86.        + 11*sin(-LS+H) -25*sin(-2*L+F) + 21*sin(-L+F);
  87.   L_MOON := P2 * frac ( L0 + DL/1296E3 ); { in rad }
  88.   B_MOON := ( 18520.0*sin(S) + N ) / ARC; { in rad }
  89.   { aequatoriale Koordinaten }
  90.   CB := cos(B_MOON);
  91.   X  := CB*cos(L_MOON);
  92.   V  := CB*sin(L_MOON);
  93.   W  := sin(B_MOON);
  94.   Y  := COSEPS*V-SINEPS*W;
  95.   Z  := SINEPS*V+COSEPS*W;
  96.   RHO:= sqrt(1.0-Z*Z);
  97.   DEC := (360.0/P2)*arctan(Z/RHO);
  98.   RA  := ( 48.0/P2)*arctan(Y/(X+RHO));
  99.   if RA<0 then RA:=RA+24.0;
  100. End;
  101.  
  102. {---------------------------------------------------------------------------}
  103. {                                                                           }
  104. { MOON: analytische Mondtheorie nach E.W.Brown (Improved Lunar Ephemeris)   }
  105. {       mit einer Genauigkeit von ca. 1"                                    }
  106. {                                                                           }
  107. {       T: Zeit in julianischen Jahrhunderten seit J2000 (Ephemeridenzeit)  }
  108. {          (T=(JD-2451545.0)/36525.0)                                       }
  109. {       LAMBDA: geozentrische ekliptikale Laenge (Aequinoktium des Datums)  }
  110. {       BETA:   geozentrische ekliptikale Breite (Aequinoktium des Datums)  }
  111. {       R:      geozentrische Entfernung (in Erdradien)                     }
  112. {                                                                           }
  113. {---------------------------------------------------------------------------}
  114. Procedure MOON;
  115.  
  116. Const PI2 = 6.283185308;
  117.       ARC = 206264.81;    { 3600*180/pi = Bogensekunden pro radian }
  118.  
  119. Var DGAM,FAC           : Real;
  120.     DLAM,N,GAM1C,SINPI : Real;
  121.     L0, L, LS, F, D ,S : Real;
  122.     DL0,DL,DLS,DF,DD,DS: Real;
  123.     CO,SI: Array[-6..6,1..4] of Real;
  124.  
  125. { berechne c=cos(a1+a2) und s=sin(a1+a2) aus den Additionstheo-  }
  126. { remen fuer c1=cos(a1), s1=sin(a1), c2=cos(a2) und s2=sin(a2)   }
  127. Procedure ADDTHE(C1,S1,C2,S2:Real;Var C,S:Real);
  128. Begin
  129.   C:=C1*C2-S1*S2;
  130.   S:=S1*C2+C1*S2;
  131. End;
  132.  
  133. { berechne sin(phi); phi in Einheiten von 1r=360 Grad            }
  134. Function SINUS(PHI:Real):Real;
  135. Begin
  136.   SINUS:=sin(PI2*frac(PHI));
  137. End;
  138.  
  139. { berechne die langperiodischen Aenderungen der mittleren Elemente }
  140. { l,l',F,D und L0 sowie dgamma                                     }
  141. Procedure LONG_PERIODIC(T:Real; Var DL0,DL,DLS,DF,DD,DGAM:Real);
  142. Var S1,S2,S3,S4,S5,S6,S7: Real;
  143. Begin
  144.   S1:=SINUS(0.19833+0.05611*T);
  145.   S2:=SINUS(0.27869+0.04508*T);
  146.   S3:=SINUS(0.16827-0.36903*T);
  147.   S4:=SINUS(0.34734-5.37261*T);
  148.   S5:=SINUS(0.10498-5.37899*T);
  149.   S6:=SINUS(0.42681-0.41855*T);
  150.   S7:=SINUS(0.14943-5.37511*T);
  151.   DL0:= 0.84*S1+0.31*S2+14.27*S3+ 7.26*S4+ 0.28*S5+0.24*S6;
  152.   DL := 2.94*S1+0.31*S2+14.27*S3+ 9.34*S4+ 1.12*S5+0.83*S6;
  153.   DLS:=-6.40*S1                                   -1.89*S6;
  154.   DF := 0.21*S1+0.31*S2+14.27*S3-88.70*S4-15.30*S5+0.24*S6-1.86*S7;
  155.   DD := DL0-DLS;
  156.   DGAM  := -3332E-9 * sinUS(0.59734-5.37261*T)
  157.             -539E-9 * sinUS(0.35498-5.37899*T)
  158.              -64E-9 * sinUS(0.39943-5.37511*T);
  159. End;
  160.  
  161. { INIT: berechne die mittleren Elemente und deren sin und cos    }
  162. {   l Anomalie des Mondes            l' Anomalie der Sonne       }
  163. {   F Abstand des Mondes vom Knoten  D  Elongation des Mondes    }
  164.  
  165. Procedure INIT;
  166. Var I,J,MAX   : Integer;
  167.     T2,ARG,FAC: Real;
  168. Begin
  169.   T2:=T*T;
  170.   DLAM :=0; DS:=0; GAM1C:=0; sinPI:=3422.7000;
  171.   LONG_PERIODIC ( T, DL0,DL,DLS,DF,DD,DGAM );
  172.   L0 := PI2*frac(0.60643382+1336.85522467*T-0.00000313*T2) + DL0/ARC;
  173.   L  := PI2*frac(0.37489701+1325.55240982*T+0.00002565*T2) + DL /ARC;
  174.   LS := PI2*frac(0.99312619+  99.99735956*T-0.00000044*T2) + DLS/ARC;
  175.   F  := PI2*frac(0.25909118+1342.22782980*T-0.00000892*T2) + DF /ARC;
  176.   D  := PI2*frac(0.82736186+1236.85308708*T-0.00000397*T2) + DD /ARC;
  177.   for I := 1 to 4 do
  178.   Begin
  179.     Case I of
  180.       1: Begin ARG:=L;  MAX:=4; FAC:=1.000002208;               End;
  181.       2: Begin ARG:=LS; MAX:=3; FAC:=0.997504612-0.002495388*T; End;
  182.       3: Begin ARG:=F;  MAX:=4; FAC:=1.000002708+139.978*DGAM;  End;
  183.       4: Begin ARG:=D;  MAX:=6; FAC:=1.0;                       End;
  184.     End;
  185.     CO[0,I]:=1.0;
  186.     CO[1,I]:=COS(ARG)*FAC;
  187.     SI[0,I]:=0.0;
  188.     SI[1,I]:=sin(ARG)*FAC;
  189.     for J := 2 to MAX do
  190.       ADDTHE(CO[J-1,I],SI[J-1,I],CO[1,I],SI[1,I],CO[J,I],SI[J,I]);
  191.     for J := 1 to MAX do
  192.     Begin
  193.       CO[-J,I]:=CO[J,I];
  194.       SI[-J,I]:=-SI[J,I];
  195.     End;
  196.   End;
  197. End;
  198.  
  199.  
  200. { TERM berechne X=cos(p*arg1+q*arg2+r*arg3+s*arg4) und   }
  201. { Y=sin(p*arg1+q*arg2+r*arg3+s*arg4)                     }
  202. Procedure TERM(P,Q,R,S:Integer;Var X,Y:Real);
  203. Var  I: Array[1..4] of Integer;  K: Integer;
  204. Begin
  205.   I[1]:=P;
  206.   I[2]:=Q;
  207.   I[3]:=R;
  208.   I[4]:=S;
  209.   X   :=1.0;
  210.   Y   :=0.0;
  211.   for K:=1 to 4 do
  212.     if (I[K]<>0) Then ADDTHE(X,Y,CO[I[K],K],SI[I[K],K],X,Y);
  213. End;
  214.  
  215. Procedure ADDSOL(COEFFL,COEFFS,COEFFG,COEFFP:Real;P,Q,R,S:Integer);
  216. Var X,Y: Real;
  217. Begin
  218.   TERM(P,Q,R,S,X,Y);
  219.   DLAM :=DLAM +COEFFL*Y;
  220.   DS   :=DS   +COEFFS*Y;
  221.   GAM1C:=GAM1C+COEFFG*X;
  222.   SINPI:=SINPI+COEFFP*X;
  223. End;
  224.  
  225.  
  226. Procedure SOLAR1;
  227. Begin
  228.   ADDSOL(   13.902,   14.06,-0.001,   0.2607,0, 0, 0, 4);
  229.   ADDSOL(    0.403,   -4.01,+0.394,   0.0023,0, 0, 0, 3);
  230.   ADDSOL( 2369.912, 2373.36,+0.601,  28.2333,0, 0, 0, 2);
  231.   ADDSOL( -125.154, -112.79,-0.725,  -0.9781,0, 0, 0, 1);
  232.   ADDSOL(    1.979,    6.98,-0.445,   0.0433,1, 0, 0, 4);
  233.   ADDSOL(  191.953,  192.72,+0.029,   3.0861,1, 0, 0, 2);
  234.   ADDSOL(   -8.466,  -13.51,+0.455,  -0.1093,1, 0, 0, 1);
  235.   ADDSOL(22639.500,22609.07,+0.079, 186.5398,1, 0, 0, 0);
  236.   ADDSOL(   18.609,    3.59,-0.094,   0.0118,1, 0, 0,-1);
  237.   ADDSOL(-4586.465,-4578.13,-0.077,  34.3117,1, 0, 0,-2);
  238.   ADDSOL(   +3.215,    5.44,+0.192,  -0.0386,1, 0, 0,-3);
  239.   ADDSOL(  -38.428,  -38.64,+0.001,   0.6008,1, 0, 0,-4);
  240.   ADDSOL(   -0.393,   -1.43,-0.092,   0.0086,1, 0, 0,-6);
  241.   ADDSOL(   -0.289,   -1.59,+0.123,  -0.0053,0, 1, 0, 4);
  242.   ADDSOL(  -24.420,  -25.10,+0.040,  -0.3000,0, 1, 0, 2);
  243.   ADDSOL(   18.023,   17.93,+0.007,   0.1494,0, 1, 0, 1);
  244.   ADDSOL( -668.146, -126.98,-1.302,  -0.3997,0, 1, 0, 0);
  245.   ADDSOL(    0.560,    0.32,-0.001,  -0.0037,0, 1, 0,-1);
  246.   ADDSOL( -165.145, -165.06,+0.054,   1.9178,0, 1, 0,-2);
  247.   ADDSOL(   -1.877,   -6.46,-0.416,   0.0339,0, 1, 0,-4);
  248.   ADDSOL(    0.213,    1.02,-0.074,   0.0054,2, 0, 0, 4);
  249.   ADDSOL(   14.387,   14.78,-0.017,   0.2833,2, 0, 0, 2);
  250.   ADDSOL(   -0.586,   -1.20,+0.054,  -0.0100,2, 0, 0, 1);
  251.   ADDSOL(  769.016,  767.96,+0.107,  10.1657,2, 0, 0, 0);
  252.   ADDSOL(   +1.750,    2.01,-0.018,   0.0155,2, 0, 0,-1);
  253.   ADDSOL( -211.656, -152.53,+5.679,  -0.3039,2, 0, 0,-2);
  254.   ADDSOL(   +1.225,    0.91,-0.030,  -0.0088,2, 0, 0,-3);
  255.   ADDSOL(  -30.773,  -34.07,-0.308,   0.3722,2, 0, 0,-4);
  256.   ADDSOL(   -0.570,   -1.40,-0.074,   0.0109,2, 0, 0,-6);
  257.   ADDSOL(   -2.921,  -11.75,+0.787,  -0.0484,1, 1, 0, 2);
  258.   ADDSOL(   +1.267,    1.52,-0.022,   0.0164,1, 1, 0, 1);
  259.   ADDSOL( -109.673, -115.18,+0.461,  -0.9490,1, 1, 0, 0);
  260.   ADDSOL( -205.962, -182.36,+2.056,  +1.4437,1, 1, 0,-2);
  261.   ADDSOL(    0.233,    0.36, 0.012,  -0.0025,1, 1, 0,-3);
  262.   ADDSOL(   -4.391,   -9.66,-0.471,   0.0673,1, 1, 0,-4);
  263. End;
  264.  
  265. Procedure SOLAR2;
  266. Begin
  267.   ADDSOL(    0.283,    1.53,-0.111,  +0.0060,1,-1, 0,+4);
  268.   ADDSOL(   14.577,   31.70,-1.540,  +0.2302,1,-1, 0, 2);
  269.   ADDSOL(  147.687,  138.76,+0.679,  +1.1528,1,-1, 0, 0);
  270.   ADDSOL(   -1.089,    0.55,+0.021,   0.0   ,1,-1, 0,-1);
  271.   ADDSOL(   28.475,   23.59,-0.443,  -0.2257,1,-1, 0,-2);
  272.   ADDSOL(   -0.276,   -0.38,-0.006,  -0.0036,1,-1, 0,-3);
  273.   ADDSOL(    0.636,    2.27,+0.146,  -0.0102,1,-1, 0,-4);
  274.   ADDSOL(   -0.189,   -1.68,+0.131,  -0.0028,0, 2, 0, 2);
  275.   ADDSOL(   -7.486,   -0.66,-0.037,  -0.0086,0, 2, 0, 0);
  276.   ADDSOL(   -8.096,  -16.35,-0.740,   0.0918,0, 2, 0,-2);
  277.   ADDSOL(   -5.741,   -0.04, 0.0  ,  -0.0009,0, 0, 2, 2);
  278.   ADDSOL(    0.255,    0.0 , 0.0  ,   0.0   ,0, 0, 2, 1);
  279.   ADDSOL( -411.608,   -0.20, 0.0  ,  -0.0124,0, 0, 2, 0);
  280.   ADDSOL(    0.584,    0.84, 0.0  ,  +0.0071,0, 0, 2,-1);
  281.   ADDSOL(  -55.173,  -52.14, 0.0  ,  -0.1052,0, 0, 2,-2);
  282.   ADDSOL(    0.254,    0.25, 0.0  ,  -0.0017,0, 0, 2,-3);
  283.   ADDSOL(   +0.025,   -1.67, 0.0  ,  +0.0031,0, 0, 2,-4);
  284.   ADDSOL(    1.060,    2.96,-0.166,   0.0243,3, 0, 0,+2);
  285.   ADDSOL(   36.124,   50.64,-1.300,   0.6215,3, 0, 0, 0);
  286.   ADDSOL(  -13.193,  -16.40,+0.258,  -0.1187,3, 0, 0,-2);
  287.   ADDSOL(   -1.187,   -0.74,+0.042,   0.0074,3, 0, 0,-4);
  288.   ADDSOL(   -0.293,   -0.31,-0.002,   0.0046,3, 0, 0,-6);
  289.   ADDSOL(   -0.290,   -1.45,+0.116,  -0.0051,2, 1, 0, 2);
  290.   ADDSOL(   -7.649,  -10.56,+0.259,  -0.1038,2, 1, 0, 0);
  291.   ADDSOL(   -8.627,   -7.59,+0.078,  -0.0192,2, 1, 0,-2);
  292.   ADDSOL(   -2.740,   -2.54,+0.022,   0.0324,2, 1, 0,-4);
  293.   ADDSOL(    1.181,    3.32,-0.212,   0.0213,2,-1, 0,+2);
  294.   ADDSOL(    9.703,   11.67,-0.151,   0.1268,2,-1, 0, 0);
  295.   ADDSOL(   -0.352,   -0.37,+0.001,  -0.0028,2,-1, 0,-1);
  296.   ADDSOL(   -2.494,   -1.17,-0.003,  -0.0017,2,-1, 0,-2);
  297.   ADDSOL(    0.360,    0.20,-0.012,  -0.0043,2,-1, 0,-4);
  298.   ADDSOL(   -1.167,   -1.25,+0.008,  -0.0106,1, 2, 0, 0);
  299.   ADDSOL(   -7.412,   -6.12,+0.117,   0.0484,1, 2, 0,-2);
  300.   ADDSOL(   -0.311,   -0.65,-0.032,   0.0044,1, 2, 0,-4);
  301.   ADDSOL(   +0.757,    1.82,-0.105,   0.0112,1,-2, 0, 2);
  302.   ADDSOL(   +2.580,    2.32,+0.027,   0.0196,1,-2, 0, 0);
  303.   ADDSOL(   +2.533,    2.40,-0.014,  -0.0212,1,-2, 0,-2);
  304.   ADDSOL(   -0.344,   -0.57,-0.025,  +0.0036,0, 3, 0,-2);
  305.   ADDSOL(   -0.992,   -0.02, 0.0  ,   0.0   ,1, 0, 2, 2);
  306.   ADDSOL(  -45.099,   -0.02, 0.0  ,  -0.0010,1, 0, 2, 0);
  307.   ADDSOL(   -0.179,   -9.52, 0.0  ,  -0.0833,1, 0, 2,-2);
  308.   ADDSOL(   -0.301,   -0.33, 0.0  ,   0.0014,1, 0, 2,-4);
  309.   ADDSOL(   -6.382,   -3.37, 0.0  ,  -0.0481,1, 0,-2, 2);
  310.   ADDSOL(   39.528,   85.13, 0.0  ,  -0.7136,1, 0,-2, 0);
  311.   ADDSOL(    9.366,    0.71, 0.0  ,  -0.0112,1, 0,-2,-2);
  312.   ADDSOL(    0.202,    0.02, 0.0  ,   0.0   ,1, 0,-2,-4);
  313. End;
  314.  
  315. Procedure SOLAR3;
  316. Begin
  317.   ADDSOL(    0.415,    0.10, 0.0  ,  0.0013,0, 1, 2, 0);
  318.   ADDSOL(   -2.152,   -2.26, 0.0  , -0.0066,0, 1, 2,-2);
  319.   ADDSOL(   -1.440,   -1.30, 0.0  , +0.0014,0, 1,-2, 2);
  320.   ADDSOL(    0.384,   -0.04, 0.0  ,  0.0   ,0, 1,-2,-2);
  321.   ADDSOL(   +1.938,   +3.60,-0.145, +0.0401,4, 0, 0, 0);
  322.   ADDSOL(   -0.952,   -1.58,+0.052, -0.0130,4, 0, 0,-2);
  323.   ADDSOL(   -0.551,   -0.94,+0.032, -0.0097,3, 1, 0, 0);
  324.   ADDSOL(   -0.482,   -0.57,+0.005, -0.0045,3, 1, 0,-2);
  325.   ADDSOL(    0.681,    0.96,-0.026,  0.0115,3,-1, 0, 0);
  326.   ADDSOL(   -0.297,   -0.27, 0.002, -0.0009,2, 2, 0,-2);
  327.   ADDSOL(    0.254,   +0.21,-0.003,  0.0   ,2,-2, 0,-2);
  328.   ADDSOL(   -0.250,   -0.22, 0.004,  0.0014,1, 3, 0,-2);
  329.   ADDSOL(   -3.996,    0.0 , 0.0  , +0.0004,2, 0, 2, 0);
  330.   ADDSOL(    0.557,   -0.75, 0.0  , -0.0090,2, 0, 2,-2);
  331.   ADDSOL(   -0.459,   -0.38, 0.0  , -0.0053,2, 0,-2, 2);
  332.   ADDSOL(   -1.298,    0.74, 0.0  , +0.0004,2, 0,-2, 0);
  333.   ADDSOL(    0.538,    1.14, 0.0  , -0.0141,2, 0,-2,-2);
  334.   ADDSOL(    0.263,    0.02, 0.0  ,  0.0   ,1, 1, 2, 0);
  335.   ADDSOL(    0.426,   +0.07, 0.0  , -0.0006,1, 1,-2,-2);
  336.   ADDSOL(   -0.304,   +0.03, 0.0  , +0.0003,1,-1, 2, 0);
  337.   ADDSOL(   -0.372,   -0.19, 0.0  , -0.0027,1,-1,-2, 2);
  338.   ADDSOL(   +0.418,    0.0 , 0.0  ,  0.0   ,0, 0, 4, 0);
  339.   ADDSOL(   -0.330,   -0.04, 0.0  ,  0.0   ,3, 0, 2, 0);
  340. End;
  341.  
  342. { Stoerungsanteil N der ekliptikalen Breite           }
  343. Procedure SOLARN(Var N:Real);
  344. Var X,Y: Real;
  345. Procedure ADDN(COEFFN:Real;P,Q,R,S:Integer);
  346. Begin
  347.   TERM(P,Q,R,S,X,Y);
  348.   N:=N+COEFFN*Y
  349. End;
  350. Begin
  351.   N := 0.0;
  352.   ADDN(-526.069, 0, 0,1,-2);
  353.   ADDN(  -3.352, 0, 0,1,-4);
  354.   ADDN( +44.297,+1, 0,1,-2);
  355.   ADDN(  -6.000,+1, 0,1,-4);
  356.   ADDN( +20.599,-1, 0,1, 0);
  357.   ADDN( -30.598,-1, 0,1,-2);
  358.   ADDN( -24.649,-2, 0,1, 0);
  359.   ADDN(  -2.000,-2, 0,1,-2);
  360.   ADDN( -22.571, 0,+1,1,-2);
  361.   ADDN( +10.985, 0,-1,1,-2);
  362. End;
  363.  
  364. { Stoerungen der ekliptikalen Laenge durch Venus und Jupiter           }
  365. Procedure PLANETARY(Var DLAM:Real);
  366. Begin
  367.   DLAM := DLAM
  368.     +0.82*sinUS(0.7736  -62.5512*T)+0.31*sinUS(0.0466 -125.1025*T)
  369.     +0.35*sinUS(0.5785  -25.1042*T)+0.66*sinUS(0.4591+1335.8075*T)
  370.     +0.64*sinUS(0.3130  -91.5680*T)+1.14*sinUS(0.1480+1331.2898*T)
  371.     +0.21*sinUS(0.5918+1056.5859*T)+0.44*sinUS(0.5784+1322.8595*T)
  372.     +0.24*sinUS(0.2275   -5.7374*T)+0.28*sinUS(0.2965   +2.6929*T)
  373.     +0.33*sinUS(0.3132   +6.3368*T);
  374. End;
  375.  
  376. Begin
  377.   INIT;
  378.   SOLAR1;
  379.   SOLAR2;
  380.   SOLAR3;
  381.   SOLARN(N);
  382.   PLANETARY(DLAM);
  383.   LAMBDA := 360.0*frac( (L0+DLAM/ARC) / PI2 );
  384.   S    := F + DS/ARC;
  385.   FAC  := 1.000002708+139.978*DGAM;
  386.   BETA := ( FAC*(18518.511+1.189+GAM1C)*sin(S)-6.24*sin(3*S)+N ) / 3600.0;
  387.   SINPI := SINPI * 0.999953253;
  388.   R     := ARC / SINPI;
  389. End;
  390.  
  391. {---------------------------------------------------------------------------}
  392. { MOONEQU: aequatoriale Mondkoordinaten                                     }
  393. {          (Rektaszension RA und Deklination DEC in Grad, R in Erdradien)   }
  394. {          T in julian.Jahrhndt. seit J2000 ( T:= (JD - 2451545.0)/36525 )  }
  395. {          Die Koord. beziehen sich auf das wahre Aequinoktium des Datums.  }
  396. {---------------------------------------------------------------------------}
  397. Procedure MOONEQU;
  398. Var L,B,X,Y,Z: Real;
  399. Begin
  400.   MOON(T,L,B,R);           { ekliptikale Moondkoordinaten           }
  401.   CART(R,B,L,X,Y,Z);       { (mittleres Aequinoktium des Datums)    }
  402.   ECLEQU(T,X,Y,Z);         { Umwandlung in aequatoriale Koordinaten }
  403.   NUTEQU(T,X,Y,Z);         { Nutation                               }
  404.   POLAR(X,Y,Z,R,DEC,RA);
  405. End;
  406.  
  407. {---------------------------------------------------------------------------}
  408. { T_FIT_MOON: Berechnet die Tschebyscheff-Entwicklung der                   }
  409. {             Koordinaten des Mondes (Reihen fuer RA,DEC und Radius).       }
  410. {                                                                           }
  411. {       TA      : Beginn des Entwicklungsintervalls (jul.Jahrh. seit J2000) }
  412. {       TB      : Ende des Entwicklungsintervalls ( TB < TA + 1 Monat )     }
  413. {       N       : Ordnung der Entwicklung                                   }
  414. {       RA_POLY,DE_POLY,R_POLY: Tschebyscheff Polynome fuer RA,DEC,R        }
  415. {---------------------------------------------------------------------------}
  416. Procedure T_FIT_MOON;
  417. Const NDIM = 27;
  418. Var   I,J,K          : Integer;
  419.       FAC,BPA,BMA,PHI: Real;
  420.       T,H,RA,DE,R    : Array[0..NDIM] of Real;
  421. Begin
  422.   if (NDIM<2*MAX_TP_DEG+1) Then WriteLn(' NDIM zu klein in T_FIT_MOON');
  423.   if (N>MAX_TP_DEG) Then WriteLn(' N zu gross in T_FIT_MOON');
  424.   RA_POLY.M := N;
  425.   DE_POLY.M := N;
  426.   R_POLY.M  := N;
  427.   RA_POLY.A := TA;
  428.   DE_POLY.A := TA;
  429.   R_POLY.A  := TA;
  430.   RA_POLY.B := TB;
  431.   DE_POLY.B := TB;
  432.   R_POLY.B  := TB;
  433.   BMA       := (TB-TA)/2.0;
  434.   BPA       := (TB+TA)/2.0;
  435.   FAC       := 2.0/(N+1);
  436.   PHI       := PI/(2*N+2);                      { h(k)=cos(pi*k/N/2)    }
  437.   H[0]      := 1.0;
  438.   H[1]      := COS(PHI);
  439.   for I:=2 to (2*N+1) do
  440.     H[I]    := 2*H[1]*H[I-1]-H[I-2];
  441.   for K:=1 to N+1 DO
  442.     T[K]    := H[2*K-1]*BMA+BPA;                { Stuetzstellen         }
  443.   for K:=1 to N+1 DO
  444.     MOONEQU(T[K],RA[K],DE[K],R[K]);
  445.   for K := 2 TO N+1 DO                          { RA stetig machen in   }
  446.     if (RA[K-1]<RA[K]) Then RA[K]:=RA[K]-360.0; { [-360,+360] !!!!      }
  447.   for J := 0 to N do                            { Tscheb.-Koeffizienten }
  448.   Begin                                         { C(j) berechnen        }
  449.     PHI  := PI*J/(2*N+2);
  450.     H[1] := COS(PHI);
  451.     for I:=2 to (2*N+1) do
  452.       H[I] := 2*H[1]*H[I-1]-H[I-2];
  453.     RA_POLY.C[J] := 0.0;
  454.     DE_POLY.C[J] := 0.0;
  455.     R_POLY.C[J]  := 0.0;
  456.     for K:=1 to N+1 do
  457.     Begin
  458.       RA_POLY.C[J] := RA_POLY.C[J] + H[2*K-1]*RA[K];
  459.       DE_POLY.C[J] := DE_POLY.C[J] + H[2*K-1]*DE[K];
  460.       R_POLY.C[J]  := R_POLY.C[J]  + H[2*K-1]*R[K];
  461.     End;
  462.     RA_POLY.C[J]:=RA_POLY.C[J]*FAC;
  463.     DE_POLY.C[J]:=DE_POLY.C[J]*FAC;
  464.     R_POLY.C[J] :=R_POLY.C[J]*FAC;
  465.   End;
  466. End;
  467.  
  468.  
  469. BEGIN
  470.  
  471. END.
  472.  
  473.